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PREFACE 


The  Army  is  interested  in  understanding  the  nature  of 
molecular  crystals  of  explosives  and  propellants.  The  long¬ 
term  objectives  of  this  research  are; 

1.  capability  to  predict  the  intrinsic  sensitivity  to  the 
initiation  of  chanical  reactions, 

2.  capability  to  predict  the  essential  physical  properties 
needed  by  new  organic  compounds  to  be  useful  as  explosives 
or  propeUants. 

There  is  a  general  feeling  that  the  advances  in  knowledge  of  microscopic 
properties  of  matter,  discovered  in  physics  and  chemistry  over  the  last 
30  years,  could  be  profitably  applied  to  this  problem,  This  research 
would  f\irnish  information  that  cannot  be  obtained  by  the  more  usual 
continuum  hydrodynamic  approach.  It  is  clear  that  several  steps  are 
involved  between  the  introduction  of  energy  into  a  solid  and  the 
eventiial  chemical  reactions  which  result.  One  area  in  which  micro¬ 
scopic  modeling  should  yield  useful  results  is  the  study  of  energy 
transfer  through  the  crystal  and  the  accumulation  of  that  energy  in  a 
particular  molecule.  It  is  obvious  that  such  steps  must  precede  the 
dissociation  of  the  molecule  and  the  formation  of  a  new  species. 

The  approach  of  our  work  is  to  bring  together  many  well  established 
ideas  from  several  branches  of  physics  and  chemistry,  making  the  minimal 
adjustments  necessary  to  use  them  on  molecular  crystals  of  explosives 
and  propellants.  Then,  we  use  these  ideas  to  describe  coherently  the 
process  whereby  ultrasonic  energy  in  a  solid  is  absorbed  into  a 
particular  molecule.  Specifically  we  look  at  the  phenomenon  of 
resonance  absorption  -  that  is,  the  process  of  energy  transfer  between 
the  lattice  vibrations  of  the  solid  and  the  lowest  internal  vibrational 
mode  of  the  molecules  which  make  up  the  solid.  We  feel  that  the 
understanding  of  chemical  reactions  caused  by  shock  or  thermal 
initiation,  as  well  as  models  for  microscopic  ’hot  spots'  caused  by 
defect^  will  be  improved  by  the  study  of  resonance  absorption. 
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displacement  coordinate  of  the  center-of-mass  from 
equilibrium  in  the  y-direction 

redefined  displacement  coordinates  in  the  y-direction  of 
the  center-of-mass  of  the  Lth  molecule  from  the  (X  ,m,n) 
equilibrium  position 

dimensionless  dvimmy  variable  defined  in  Equation  (35d) 

displacement  coordinate  of  the  center-of-mass  from 
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pmoportional  to  the  second  derivative  of  U(X,Y,Z)  evaluated 
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frequency  of  the  sound  wave  in  Table  6 
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maximum  lattice  frequency 

the  Lennard-Jones  parameter  for  the  distance  of  closest 
approach  of  two  molecules  which  collide  with  zero  initial 
relative  kinetic  energy 


T 

l’i> 


0) 


“d 


relaxation  time  equal  to  l/f^^ 

defined  in  Equation  (10c) 

eigenfunction  of  defined  in  Equation  (23) 

circular  frequency  for  the  lowest  vibrational  mode  of  the 
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Debye  circular  frequency  defined  in  Equation  (28b) 

circular  frequency  for  the  jth  branch  of  the  mode  with 
wavevector  k,  defined  in  Equation  (14) 

maximum  circular  lattice  frequency  defined  in  Equation  (15) 
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I.  INTRODUCTION 


Originally,  we  became  interested  in  studying  the  way  that  lattice 
vibrational  (phonon)  energy  is  absorbed  into  the  molecular-vibrational 
node  of  a  defect  molecule  or  a  molecule  located  next  to  a  defect.  In 
these  organic  crystals  the  forces  between  the  atoms  in  a  molecule  are 
much  stronger  than  the  van  der  Waals  -  type  forces  between  molecules 
in  the  crystal.  The  energy  necessary  to  excite  hhe  0->-l  transition  in 
the  lowest-frequency  mode  of  a  molecular  vibration  is  usually  much 
higher  than  the  energy  of  the  highest-frequency  phonon  in  a  perfect 
crystal.  Thus,  a  multiphonon  process  is  usually  required  to  excite 
a  molecular  vibration.  During  the  course  of  our  literature  search  for 
our  previous  paper!,  we  found  some  papers  by  Liebermann^"^  which 
coupled  the  lattice  vibrations  to  the  molecular  vibrations  in  perfect 
crystals  in  essentially  a  standard  way  and  suggested  that  the  ultra¬ 
sonic  absorption  observed  in  benzene  could  be  explained  by  this 
coupling.  Other  investigators^-S  have  subsequently  observed  ultrasonic 
absorption  in  other  crystals  and  used  Liebermann's  model  to  analyze 
their  results.  As  the  agreement  between  theory  and  experiment  left 
some  doubt  as  to  whether  or  not  Liebermann's  model  adequately  explained 
the  phenomenon  called  'resonance  absorption*,  we  began  to  inspect  his 
model  very  carefully.  We  have  discovered  that  Liebermann  made  an 
unfortunate  error  in  carrying  out  the  details  of  his  model,  which  in¬ 
validates  his  final  theoretical  expression  for  the  transition  rate  for 
lattice  energy  to  be  absorbed  into  the  molecule.  As  this  expression 
is  still  being  used  by  the  other  investigators,  we  have  decided  to 
formulate  a  model,  similar  to  the  one  of  Liebermann,  and  calculate  an 
expression  for  the  transition  rate.  Our  new  result  can  be  assessed  by 
the  previous  investigators  on  their  experimental  work,  and  the  direction 
of  future  research  on  this  phenomenon  can  be  suggested  from  such  studies. 

In  our  model  we  consider  a  simple  cubic  lattice  composed  of 
molecules  separated  by  a  constant  lattice  spacing.  We  make  the  approxi¬ 
mation  that  the  center-of-mass  motion  of  each  molecule  is  nearly  inde¬ 
pendent  of  the  motion  of  the  atoms  with  respect  to  the  center  of  mass. 


1.  D.F.  Strenzwilk  and  P.D.  Yedinak,  Phys.  Rev.  7.>  874  (1973). 

2.  L.  Liebermann,  Phys.  Rev.  115,  1052  (1959). 

3.  L.  Liebermann,  J.  Acoust.  Soc,  Am.  1073  (1959), 

4.  L.  Liebermann," Physical  Acoustics",  Vol,  IV,  Part  A  (Acad,  Press, NY, 
1966)  ed.  by  W.P.  Mason,  p.  183. 

5.  R.A.  Rasmussen,  J.  Chem.  Phys.  1821  (1962). 

6.  S.S.  Yun,  and  R.T,  Beyer,  J.  Chem.  Phys.  4^,  2538  (1964). 

7.  J.D.  Wilson,  and  S.S.  Yun,  J.  Acoust.  Soc.  Am.  164  (1971). 

8.  R.A.  Rasmussen,  J.  Chem.  Phys.  4^,  211  (1967), 

9.  A.E.  Victor,  and  R.T.  Beyer,  J.  Acoust.  Soc.  Am,  54,  1639  (1973). 
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Then  the  term  which  couples  the  two  motions  contains  the  only  inhar¬ 
monic  factors  in  our  model .  The  center-of-mass  motion  is  treated  in 
the  same  manner  as  a  monatomic  simple  cubic  lattice  in  which  the 
nearest  neighbors  interact  through  a  central  force  derivable  from  a 
Lennard-Jones  '•6-12"  potential  and  a  non-central  force  which  we  have 
arbitrarily  made  consistent  with  the  Debye  frequency.  The  molecular 
vibrations  are  assumed  to  be  the  same  as  in  the  gas  phase.  Using  a 
modified  Debye  distribution^ifor  the  lattice  frequencies  and  the  Planck 
distribution  for  the  occupation  numbers  of  the  various  modes,  we  have 
calculated  the  transition  rate  for  a  molecule  to  absorb  energy  from  the 
lattice. 


II.  THEORY 

A.  Coupling  Term. 

Consider  the  Lth  molecule  in  a  simple  cubic  crystal  with  its  six 
nearest  neighbors  at  equilibrium  lattice  separation,  d,  and  assume  a 
Lennard-Jones  "6-12"  potential  between  molecules.  If  the  center  of 

mass  of^Lth  molecule,  whose  equilibrium  position  is  denoted  by  = 

d(Ci  +  mj  +  nk),  is  moved  a  small  distance  (X,Y,Z),  while  the  nearest 
neighbor  molecules  are  held  fixed,  (as  the  other  investigators2-8  have 
done  for  simplicity),  we  can  write  the  potential  energy  for  the  Lth 
molecule  in  terms  of  the  molecule  located  at  (0,0,0)  without  loss  of 
generality  as  follows. 


CD 


where  n  and  a  are  constants  depending  of  the  kind  of  molecules  involved. 
2  2  2 

and  C  +  ra  +  n  =  1,  and  E ,m,n  =  0,1,-1. 


If  we  expand  this  potential  in  a  Taylor  series  to  the  fourth 
order  about  the  origin  (0,0,0)  in  X,Y,Z,  we  obtain: 

UCX,Y,Z)  =  U(0,0,0)  +  Xj(X^  +  Y^  +  Z^)  +  +  y"^  +  z"^) 

+  x^cxV  +  xV  +  Y^Z^)  + 


(2) 
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where 
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also  have  to  take  into  account  the  forces  between  atoms  in  different 
molecules  and  not  just  the  average  intermolecular  force  between  centers. 
An  example  of  a  more  rigorous  treatment  of  this  coupling  can  be  found 
in  Rasmussen's^  earlier  paper  on  solid  cyclohexane,  but  in  a  later 
paper®  on  carbon  tetrachloride  he  used  the  simpler  form  of  Liebermann. 

We  also  assiune  that  the  amplitude  of  the  center  of  mass  is  much  greater 
than  the  molecular-vibrational  coordinate,  i.e.. 


^mn 

^£mn 

^£mn 

» 

^£mn 

max. 

max, 

max 

If  we  substitute  Eq.  (4)  in  Eq.  (2),  neglecting  not  only  terras  which  do 
not  mix  the  coordinates  x,y,z,  and  v  since  they  contribute  to  the  un¬ 
coupled  energy,  but  also  terms  involving  x^v^,,..,  and  xv®,...  because 
of  Eq.  (5),  we  have  the  coupling  between  the  lattice  and  molecular 
vibrations,  viz.. 


2  4  3  3  3 

^int  ~  ^  ^l^£mn^^mn*^£mn*^fimn^  *  ^  ^2^£mn^^mn*^£mn*^£mn^ 


2  2  2  2 
mn^£  mn^^£  mn^  mn^^  mn^£  mn^^£  mn^  mn 


(6) 


*  ^fimn^£iin^^fimn^£mn^  * 

where  the  interaction  Hamiltonian  or  coupling  term  is 
B.  Total  Energy  of  Crystal 

The  main  contribution  to  the  energy  of  this  system  is  composed 
of  the  lattice  dynamical  energy  of  the  molecules  as  a  whole,  and  the 
energy  of  the  molecular  vibrations.  For  the  lattice  energy  in  the 
harmonic  approximation,  we  have  chosen  the  model  of  Montroll -Potts-Rosen- 
stock-Newell  (MPRN)^0»ll  where  the  lattice  coordinates  (x,y,z)  appear 
quadratically  in  the  potential  energy.  The  Hamiltonian  for  the  lattice. 


^®E.W.  Montroll,  and  R.B.  Potts,  Phys.  Rev.  100,  525  (1955). 

^^H.B.  Rosenstock  and  G.F.  Newell,  J.  Chem,  Phys.  21,  1607  (1953). 


18 


^lat  “  ^  ^  sin  ^^inn''’^8inn*^8inn^  *  ^  \  £ran^  ^*8inn  ^+l,Tnn^ 

*  ,m+l,n^  *  n+1^  ^ 

*  ^  ^nc  £m  ^  ,m+l,n^ 

2  2  2 

*  ^^«Bm"^Ein,n+l^  *  ^^Emn'^+l.mn^  *  ^^]rai"^m,n+l^ 

< 

where  M  is  the  mass  of  a  molecule,  =  2Xj^  is  the  central  force 

constant  and  y  is  the  noncentral  force  constant, 
nc 

For  a  perfect  lattice  composed  of  identical  molecules  of  mass  M, 
the  equations  of  motion  fora  system  composed  ofNxNxN  =  N 
particles  with  nearest -neighbor  central  and  noncentral  interactions 
only  are 

^  ^  mn  “  \  ^^+l,mn‘^^mn'^^-l,mn^  *  ^nc^^  ,m+l,n‘^^mn'^^  ,m-l,n^ 

*  ^nc^^m,n+l'^^mn^^m,n-l^  *  ^  ^ 


^  ^  mn  \^£  ,m+l,n"^^£mn‘^^fi,m-l,n^  *  ^nc^^£m,n+l‘^^£mn‘^^fim,n-l^ 


*  \cfy£  *l,mn-2y£mn^y£-l.mn^  ' 


(8b) 


and 


•  • 

^  ^  mn  "'^c^^£m,n+l  ^^£inn'*'^fim,n-l^  ^nc^^fi+l,mn  ^^£mn  £-l,mn^ 


*  ^nc^*£,m+l,n  ^^fimn*^£  ,m-l,n^  ' 

£,m,n  =  -jN+l,  ...  ,yN,  N^ven  ’ 
where  we  assume  periodic  boundary  conditions 
*£mn  ~  ^+N,mn  ~  ^  ,m+N,n  ~  ^m,n+N  ’ 


(8c) 


(9) 
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The  normal  modes  for  Eq.  (8a)  can  be  expressed  in  terms  of  the  ampli¬ 
tude  A  as 

X 

^mn  ”  \  exp(-ia)t)  ,  (10a) 

where 

u02mn)  =  exp  i((|)jfi  +  (|)2  m  +  (|)j  n)  (10b) 

and  for  j  =  1,2,3 

<()j  a  (2Traj/N);  -  1,2,...,N  (10c) 

or  -  J2N+l,...,JiN. 

Similar  solutions  satisfy  Eqs.  (8b)  and  (8c). 

12  ... 

The  general  solution  is  given  by  the  following  expansion  in 

terms  of  plane  waves: 


mn 


(11) 


j=l,2,3 


where 


•  Q  (j)  i 

•  ^  _  n _  ^1-  ^  _ I- 


is  the  complex, mass-weighted,  normal 


:  -..(y 


coordinate  for  the  jth  branch  of  the  mode  with  wavevector  R 
is  the  x-componeiit  of  the  vector  e  with  the  definition  in  the 
MPRN  model  '  / 

i  =  (1,0,0);  i  ^2)  =  (0»1.0);  e  =  (0,0,1)  ,  (12) 

and  wavevector  k  is 


ic 


2ir 

Nd 


1  +  a2  j 


+  a3  k)  . 


(13) 


^^A.A.  Maradudin,  E.W.  Montroll,  and  G.H.  Weiss,  "Solid  State  Physics", 
supp.  3,  (Acad.  Press,  NY,  1963). 
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The  y  and  t  displacements  have  similar  expansions.  Corresponding  to 


each  mode^^^} is  a  circular  frequency  which  satisfies 
,  3 

M  wTC  £  )  =  2  Y,  Cl  -  cos  <().)  , 

J  i=l  ^  ^ 


where  y.  =  y  and  y-  ®  Y 
'j  'c  nc 


The  maximum  lattice  frequency  is 


2.,  V 

0)  »  4  - 

°  M 


C.  Quantum  Mechanics 


The  energy  of  the  lattice  in  Eq,  (7)  can  be  quantized  by 
expressing  Q  j  by  the  following  transformation  which  can  be  found 
'  '  12 

in  standard  texts  : 


'  \  2a),  (k)  / 


Cl  6a) 


'^(i)  ”’0)  “ 


where  a  and  are  the  creation  and  annihilation  operators, 
respectively,  for*' phonons  whose  wavevector  is  ic  and  which  belong  to 
the  jth  branch,  and  h  is  Planck’s  constant  divided  by  2it.  These 
operators  obey  the  usual  commutation  relations  of  lattice  dynamics 
and  have  familiar  matrix  elements  for  Bose  operators  in  the  numbers 
representation.  Thus,  1^.  (7)  becomes  a  sum  of  quantum-mechanical 
operators  in  the  numbers  representation  when  Eq,  (16)  is  substi¬ 
tuted  into  either  Eq.  (11)  or  the  similar  expressions  for  the 
y  and  z  displacements,  i.e.. 
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(17) 


H,  =  S  h  0) .  (ic)  (a*,-  .Oir  . 
lat  3  k3  kj 

j 


+  1/2)  . 


The  molecular-vibrational  energy  for  the  lowest  vibrational  mode 
of  the  Lth  molecule  in  the  crystal  can  be  expressed  in  terms  of  its 
reduced  mass  the  vibrational  frequency  Q,  and  the  normal  coordinate 

for  this  mode  v(L),  so  that  the  total  vibrational  energy  can  be  expressed 
as 


L 


(18) 


where  H  ,  is  the  molecular-vibrational  Hamiltonian.  This  Hamiltonian 
mol 

can  be  quantized  by  expressing  the  normal  coordinate  v(L)  by  the 
following  trans format ion? 


and 


(19a) 


(19b) 


where  b^  and  b.  are  the  creation  and  annihilation  operators,  respectively, 
L  ij 

These  operators  obey  the  usual  commutation  relations  and  have  the 
familiar  matrix  elements  for  Bose  operators  in  the  numbers  representa¬ 
tion.  Thus,  Eq,  (18)  also  becomes  a  sura  of  harmonic  oscillator 
operators,  one  for  each  molecule,  viz.. 


H  =  S  hQ  (b.^  b,  +  1/2)  . 
mol  "  '•  L  L  ' 


(20) 


The  total  uncoupled  energy  of  the  crystal  in  our  model  is  given  by 
the  addition  of  the  lattice  energy  of  Eq.  (17)  and  the  molecular- 
vibrational  energy  of  Eq.  (20)  to  yield  the  Hamiltonian, 

%  =  "lat  *  "mol  • 

which  operates  on  the  eigenket  \^>  that  is  the  product  of  harmonic 

oscillator  functions  in  the  niMbers  representation,  i,e,. 
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-  Eq1\|;>  , 


where 


lt|>  -  ln(l).n(2),....n(N^);  n(j^)  ,  n(3^),...,n(^^)>,  (23) 


where  n(l),  n(2),  etc.  are  the  molecular-vibrational  quantum  numbers  cor‘* 


responding  to  the  Lth  molecule,  and  the 


"(3 


are  the  phonon  occupation 


numbers  corresponding  to  the  lattice  vibrational  mode  designated  by  the 
wave  vector  R.  and  acoustical  branch  j.  There  are  no  lattice  optical 
branches  in  t^e  MPRN  model. 

. . '  _ _ _  3 

The  eigenvalue  Eq  is  the  sum  of  4N  harmonic  oscillators,  i.e., 

Ert  =  E  (n(  JV  1/2)  h  Oj^ic)  +  E  Cn(L)  +  1/2)  h  .  (24) 

^  k,j  VJ  J  L 

If  we  quantize  Eq.  (6)  by  the  same  procedures  used  for 
Eq.  (17)  and  Eq.  (20) ,  becomes 


H.  ^  =  —  A 


*  5  **  0) 


(2mn) 


/k\ 

where  Q  K  I  and  v(L)  are  the  operators  defined  in  Eq.  (16)  and  Eq.  (19), 

respectively.  Hereafter,  the  first  term  in  Eq,  (25), which  can  transfer 
energy  only  when  the  phonon  frequency  equals  the  molecular  vibrational 
frequency,  will  be  dropped  since  £2>a)^  in  our  model. 

The  total  Hamiltonian  operator  H^^^  is  the  sum  of  Eq,  (21)  and 
Eq.  (25),  i.e. 

tot  0  int* 


(26) 


where  contains  the  only  anharmonic  terms  in  this  model.  It  is 
this  interaction  Hamiltonian  H.  which  allows  the  energy  to  flow  from 

^11  t 

the  lattice  into  a  molecule.  From  time-dependent  perturbation  theory 
the  probability  of  going  from  the  initial  state  |  i>  to  the  final 

state  |.£>,  when  the  lattice  and  tth  molecule  at  time  t  =  0  are  in  a 
state  that  exactly  corresponds  to  one  of  the  eigenstates  of  the 
unperturbed  system  |  i>,  is 


where  the  Dirac  6  function  is  used  and  the  superscripts  f  and  i  refer 
to  the  final  and  initial  states,  respectively. 


There  are  many  final  states  | £>  in  which  the  system  could  be  after 
making  a  transition  from  the  given  initial  state  |i>.  The  density  of 
states  per  frequency  g((Oj(k))  has  been  worked  out  for  this  model  by 

Montroll^*^*^^.  He  showed  that  the  Debye  frequency  distribution  was 
given  as 


g(a).(k))  =  3N^  a)j(R)/a)p 


0  <  u  j  (k)  <  ojp 


=  0 


(0^  (ic)  >  tOp 


(28a) 


where  the  circular  Debye  frequency  Up,  (which  is  2^  times  the  Debye 
frequency  v^,  i.e.,  u)p=2TfVp),  equals 


0) 


3 

D 


(28b) 


We  decided  to  renormalize  the  Debye  distribution  so  that  the  maximum 
frequency  a)^=2Trv^  of  the  MPRN  model,  which  is  greater  than  oOp,  would 

be  included  in  our  frequency  spectrum,  thus 


g(a).(k))=  3N^  o)^(k)/a)^ 
=  0 


0  <  0) .  (ic)  <  0)^ 

3  0 

Wj  (t)  >  Wq 


(29) 


13 

E.W.  Montroll,  in  "Third  Berkeley  Symposium  on  Mathematical  Statistics 
and  Probability",  ed.  by  J.  Neyman,  (Univ.  Cal.  Press,  Berkeley  and 
Los  Angeles,  1956),  p.  209. 
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The  purpose  of  the  renormalization  was  that  we  wanted  to  keep  the 
high  frequency  phonons  in  our  distribution  which  are  truly  present, 
but  we  did  not  want  to  add  the  additional  mathematical  difficulty  of 
using  the  exaet  density  of  states  which  can  be  found  in  Montroll's  paper. 
The  transition  rate  R  gives  the  total  transition  probability  per  unit 
time  that  the  system  initially  in  the  state  |  i>  is  in  any  of  the  many 
final  states  I  f>,  i.e.. 


Pjf  8(»j(£i))  gt«j,tl2))  gCiSj/zClj)) 


do)j(kj)  dcOj,  (^2)  du)j//Cic^) 


where  R  is  dependent  on  a  three-phonon  process  involving  the  lattice 
mode  frequencies  ),  a).,/(k„)  and  the  sum  is  taken  over 

the  nine  allowed  values  of  the  j's  in  Eq.  (25). 


lantum  Statistics 


In  Eq.  (23)  we  wrote  the  eigenket 


occupation  numbers 


|t()>  in  terras  of  phonon 


As  far  as  quantum  mechanics  is  concerned 


these  numbers  are  integers  and  can  be  selected  arbitrarily,  but  for  a 
system  in  thermal  equilibrium  at  temperature  T  we  know  that  it  is 
harder  to  populate  the  higher  energy  levels  according  to  quantum 
statistics.  Therefore,  for  the  phonon  occupation  numbers  for  the 
initial  state  li>  we  have  used  the  Planck  distribution  where  the 


l-C')! 


*s  are  chosen  according  to 


n(j^)>  =  ^exp(ha).(kp/kg  T)  -  1  ^ 


where  <. .  .>  denotes  the  average  in  thermal  equilibrium  and  kg  is  the 
Boltzmann  constant. 


III.  CALCULATION 


Transition  Rate 


To  calculate  the  transition  rate  R  of  Eq.  (30)  we  must  calculate 

the  matrix  elements  !<  flH.  li  >1^  of  H.  of  Eq.  (25)  that  must  be 

'  '  int '  '  int 

inserted  into  of  Eq.  (27).  Let  us  first  consider  the  stimmation 
over  the  allowed  values  of  the  j's  in  Eq.  (30).  For  a  given  kj,k2,R2 
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there  are  nine  non-zero  matrix  elements  describing  the  transition  from 
the  initial  state  |i>  =  InTLlrnl  .*1  .  nl  .71  .  nl  .7»l  >to  final 


states 


>0  matrix  elements  aescrioing  rne 

-  |n(L,;n(Jl)  ,  „('?)  .  n(;?,) 

|£>  -  |nCL)  •*  1;  -  1.  -  1,  "(j//)-  »• 


We  have 


adopted  the  convention  of  writing  the  occupation  numbers  of  only  those 
states  that  will  change  during  the  transition.  Thus,  |i>  is  really  of 
the  form  |i|;>  in  Eq.  (23),  where  the  most  probable  occupation  numbers 
are  given  by  Eq.  (31).  The  matrix  elements  for  the  following  cases 


are: 


Case 


1.  j=  y  *  j";  If  >  =  |nCL)+l;  -1,  -1,  n(j^) 


|<  flH.  Jd>|2  =  4  (  *  3/'2  )  (  - - ^ 


1> 


(„(LM)  n ('*)„( '2)  n{})  . 


(32) 


There  are  three  terms  of  the  above  form  for  j=l,2,3. 

Case  2.  j  j'  =  j'' ;  If  >  =  |n(L)+l; 

1<  flH.  Ji>l2  =  4  (  ^  3  5V2  )  ( - - - ~r -  ) 

4^  \  /  \  y  ojjCtj)  a)j/Ck2^  ^  / 


(n(L) 


(33) 


There  are  six  terms  of  the  above  form.  When  we  substitute  these  nine 
terms  into  Pj^£  appearing  in  R,  Eq.  (30) ,  along  with  the  proper 

g(a).(ic))  of  Eq.  (29),  we  find  that  there  are  3  similar  integrals  for 

Case  1  and  6  similar  integrals  for  Case  2.  Therefore,  the  form  of  R 
in  Eq.  (30)  becomes 
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(34) 


^  ^  a  xyz^^  " ^ ^ 

where  the  integral  has  been  made  dimensionless  by  using  the  following 
definitions 


c.ii  1 

^  m3 
VM 


a  =  J2/ujq  , 


(1/2  +  1/4  X^  ) 

(2Tr)^ 


(?)  ■ 


(35a) 


a  -iiWQ/kgT  =  e/T  , 


(35b) 
(35c) 

u)(kj^)  =  X  (Up,  u)(k2)  =  y  oJq,  (oCkp  =  z  oip,  (35d) 

n(kj)  =  n(ax)  =  (exp(ax)  -  1)'^,  nCk^)  =  nCay),  nCk^)  =  n(az),  (35e) 
and 

2  .  .  2  9 

(35f) 


g(x)  =  x^,  g(y)  =  y^,  g(z)  =  z^  . 


The  domain  of  integration  of  Eq.  (34)  must  be  determined  subject  to  the 
condition  that 

0  <  x,y,z  <  1,  (36) 

since_the  lattice  fre<^uencies  must  lie  between  0  and  to  or  else 
g(a)j(k))  in  Eq.  (29)  will  be  zero.  ° 

If  we  integrate  over  the  variable  z,  the  delta  function  restricts 
us  to  the  condition 


z  =  a-x  -y  . 

Eq.  (36)  gives  the  following  condition 
0  <  a-x  -y  <  1, 
from  which  we  find  that 

y  >  -X  +  (a  -  1),  (lower  limit) 

and 

y  <  -X  +  a.  (upper  limit) 


(37) 

(38) 

(39) 

(40) 
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Therefore,  in  the  xy  plane  we  must  integrate  over  all  values  of  (x,y) 
that  fall  between  the  upper  and  lower  limits  and  are  common  to  the 
unit  square  in  the  first  quadrant.  In  our  model  the  ratio  of  the 
molecular  vibrational  frequency  to  the  maximum  lattice  frequency  lies 
between  2  and  3  for  benzene,  p-didilorobenzene,  carbon  tetrachloride,  and 
cyclohexane,  i.e.,  2  <  a  <  3,  whereas  for  naphthalene  the  ratio  lies 
between  1  and  2,  i.e.,  1  <  a  <  2. 

The  domains  of  integration  for  these  two  cases  can  be  seen  in 
Figure  1.  However,  the  formula  for  R  developed  by  Liebermann^”^ 
integrates  over  the  whole  unit  square  without  restricting  the  domain 
of  integration.  According  to  our  derivation  such  an  integration  would 
not  conserve  energy.  In  the  model  of  Liebermann 

a  =  —  =  3.90, 

% 

since  the  lowest  molecular  vibrational  frequency  of  benzene  is 

a  =  2irfl.21  X  lO^^s"^)  and  the  maximum  lattice  frequency  was  arbitrarily 

^  12  **1 
chosen  to  be  the  Debye  frequency,  ^  2tt(3.1  X  10  s  In  that 

case  the  lower  limit  for  the  integration  domain  is  the  line  y  =  2,9-x, 
which  passes  entirely  above  the  unit  square.  The  integration  domain 
is  therefore  empty  and  consequently  R  =  0,  Liebermann *s  formula 
yields  a  non-zero  value  of  R  for  a  -  3,90  and  must  be  in  error. 

There  is  no  way  that  three  phonons,  whose  maximum  possible  frequency 
is  (Dq,  can  give  the  energy  necessary  to  excite  a  molecular  vibrational 
mode  of  frequency  ^  =  3.90  cdq. 

The  explicit  form  for  R  of  Eq.  (34)  for  the  case  2  <  a  <  3  is 


R  =  C  I(a,a), 

aa 


(41) 


where 


I  (a 


-  ij 


(42) 


g'^  xy(a.-x-y)  dy  dx 


x=a.-2  y=-x+Ca-l) 


Ce”*  -  11  Ce'‘>'  -  1)  -  1) 


We  decided  that  an  analytical  expression  for  the  function  of  Eq, 
(42)  would  be  of  more  general  use  than  a  niomerical  solution  on  the 
computer  for  a  few  specific  values  of  (t  and  a.  In  order  to  arrive  at 
such  an  expression,  we  used  the  Taylor  series  expansion. 


ax 


+  B. 


(oxii 

2! 


r  i4 

(ax) 

4! 


+  B. 


(ax)^ 

6.' 


(ax)^ 

8! 


+ 


•  ,  ,  , 


|axl<  Zir  , 


(43) 


y= “X+^ 


yS“X+(^-l) 


2<a<2 


[\,a-2) 


where  the  B  's  are  the  Bernoulli  numbers,  (B,  *  1/6,  B„  =  1/30, 
n  1  z 

Bj  =  1/42,  B^  =  1/30,  etc.)*  Whbn  we  substituted  Eq.  (43)  and  the 

similar  expressions  for  ay  and  a(a-x-y)  into  Eq.  (42),  we  obtained 
a  power  series  in  a  of  the  form 


I  (a,  a) 


00 


n=0 


0^(0.)  a  , 


(44) 


where  we  evaluated  c  for  n  =  0  to  n  =  9.  The  first  few  c  's  were 

n  n 

calculated  by  hand,  but  as  this  task  became  exceedingly  tedious  we 

relied  much  more  on  the  values  generated  by  a  symbol  manipulation 

program  developed  by  George  C.  Francis^^  of  BRL.  In  the  Appendix  A 

of  this  BRL  Report,  we  show  the  integral  expressions  for  the  c^'s 

which  will  be  given  here  in  terms  of  a.  only  in  Table  1. 


It  is  at  this  part  in  our  paper  that  we  want  to  point  out  our 
disagreement  with  another  approximation  which  Liebermann  made.  He 
truncated  the  series  expansion  of  Eq.  (43)  after  the  first  term.  But 
in  our  calculation  for  benzene,  (as  we  will  show  later),  a  =  .965  in 
Eq.  (35)  at  a  temperature  of  273K,  and  the  numerical  value  of  I  (cl, a) 
of  Eq.  (44)  changes  appreciately  by  retaining  higher  powers  of  a. 


An  alternate  expression  for  I  (a, a)  in  Eq.  (42)  was  pointed  out 
to  us  by  Walter  Egerland^^  of  BRL.  By  making  three  successive 
substitutions,  the  variable  lower  limits  of  integration  are  replaced 
by  zero.  This  form  of  I(a,o)  is  better  suited  to  numerical  integration 
on  a  computer,  and  also  shows  how  the  factor  (3-a)2  in  Eq.  (44) 
naturally  occurs.  First,  substitute  y-^t-x+a-1  in  the  y  integration 
and  then  substitute  t-»-(x-a+2)v.  The  resulting  integration  goes  from 
v=0  to  1.  Finally,  substitute  x->-(a-2)  +  (3-a)  u  in  the  x  integration. 

The  resulting  integration  goes  from  u»0  to  1.  Therefore,  I(a.,o) 
becomes 


I  (a, a)  = 


3  2 

acs-ay 


f  lilEl 

L 


p+(l-p)u] 


0 


(45) 


f  [  (l-p)uv-(l-p)u+ll  r  l-(l-p)uv]dv  du 

^  4=0  (e«t  (l-P)«v-(l-p)u+l]_^)  (^a[l-(l-p)uv]_j)  ’ 

^^George  C.  Francis,  "SYMAP2  -  An  Operational  Computer  Based  Algebraic 
Symbol  Manipulator'*,  BRL  Memorandtim  Report  No.  2199,  (AD  746  957), 
June  1972. 

^^Walter  Egerland,  private  communication. 
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TABLE  1.  VALUES  OF  THE  Cn's  ARE  LISTED  AS  A 

FUNCTION  OF  a=  a/u)^. 

Co  1/2 
Cl  -a/A 

C2  (5a^+2a-3)/96 
C3  -(a^+2a^-3a)/192 
C4  (13a^+162a^-165a?-228a+162)/86,400 
C5  (a^-6a^-30a^+114a2-81a)/86,400 

C*  (67a*-2910Al3,941a'*-21J20a.^-9.324a?  +  40,140a-17,838)/152,409,600 

C7  (-23a^+132a*+747a^-6,198a^+16,569a^-20,070a^+8,919a)/152,409,600 

Cg  (-364A8706d^-60,318a*+184,218<i^ -234, 270a^-58,608a^+533,637a?-562, 266a +  176, 931) 
/2.286144  X  10’° 

C9  (lla’-74a®-448«^+6058a*-26,610a^+62,262a^-84,068a^+62,474a2-19,659a) 

/5.08032  X  10’ 
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where 


p  ■  a  -  2,  u  =  (x-p)/Cl-p),  and  v  =  (y+x-p-l)/(x-p) ,  (46) 


An  expansion  of  I (a, a)  in  power  series  of  a-S  is  shown  in  Appendix  B. 

For  values  of  a  <  1  we  obtained  numerical  agreement  with  Eq,  (44)  by 
using  the  first  four  non  vanishing  terms  of  that  expansion.  However, 
we  prefer  Eq.  (44)  where  the  temperature  dependence  is  in  the 
foirm  of  a  power  series  in  T"^  rather  than  complicated  exponential 
factors  for  the  latter  case. 

B.  Application  to  Real  Crystals 

Up  to  this  part  of  our  paper  we  have  described  how  to  calculate 
the  transition  rate  for  a  hypothetical  molecular  crystal  with  simple 
cubic  structure  and  central  and  noncentral  forces.  Now  we  will 
discuss  the  approximations  we  used  to  treat  actual  molecular  crystals. 
Liebermann2“4  assumed  that  the  crystal  structure  of  benzene,  which  is 
simple  orthorhombic,  could  be  approximated  by  a  simple  cubic  structure 
with  the  length  of  the  cube,  d,  assumed  to  be  the  nearest  neighbor 
distance.  We  have  made  the  same  assumption  for  benzene,  as  well  as 
for  cyclohexane^  and  carbon  tetrachloride®,  which  are  face-centered 
cubic,  and  p-dichlorobenzene®  and  naphthalene®*',  which  are  monoclinic 

C.  We  have  used  the  data  contained  in  the  above  papers  and  listed  in 
Table  2  to  calculate  the  central  force  constant  Y^,  i.e,, 

Y  =  2X,  (47) 

c  1 

where  is  defined  in  Eq.  (3a).  The  results  are  given  in  Table  3, 

For  the  noncentral  force  constant  we  were  faced  with  the  follow¬ 
ing  dilemma.  We  had  no  information  on  what  noncentral  forces  to  use, 
and  yet  if  we  used  only  central  forces  and  a  simple  cubic  lattice,  we 
would  have  three  independent  linear  chains  in  the  x,y,  and  z  directions. 
Our  model  would  reduce  to  a  one  dimensional  problem  with  the  density 
of  states  g(ujj(k))  the  same  as  for  the  linear  chain.  Therefore,  we 

decided  to  determine  the  noncentral  force  constant  Yj^^.  arbitrarily 

by  making  it  consistent  with  the  known  Debye  frequency  Vp  as  reported 

in  the  aforementioned  papers.  Using  Eq.  (28)  with  the  relation 

“d  “  relation  we  determined  Y^^^. 

and  the  maximum  lattice  frequency  to  be  used  in  Eq.  (29).  The  results 
are  listed  in  Table  3.  To  determine  X^^,  X2,  and  X^  of  Eq.  (3),  we 

used  the  information  in  Table  2  and  listed  the  results  in  Table  4, 
Finally,  using  values  of  k_  =  1.3806  X  10”l®  erg-K"l  and  h  -  1,0542 
X  10“^'  erg-s,  we  calculated  values  of  the  constant  C  in  Eq.  (35a)  for 
the  various  substances  and  also  evaluated  the  constant  fl.  in  Eq.  (35b) 


32 


Table  2.  Constants  Compiled  From  Other  Investigator's  Papers. 
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and  the  ratio 


=  0  in  Table  5. 


We  now  have  enough  information  to  check  our  results  for  benzene 
p-dichlorobenzene^,  carbon  tetrachloride®,  and  cyclohexane^  against 
those  of  previous  investigators,  who  measured  the  absorption  at  the 
temperatures  indicated  in  Table  6,  and  determined  Rg™^from  a  gas  phase 
expression  used  by  Liebermann,  The  expression  for  the  absorption  of 
ultrasound  given  by  Liebermann2~4  is 


0)^  R 


'•  expt 


where  0  is  the  absorption,  V  the  velocity  of  sound,  C^,  are  the 

specific  heats  at  constant  pressure  and  constant  volume, is  the 
specific  heat  associated  with  that  part  of  the  molecular  energy  which 
relaxes,  u  is  the  circular  frequency  of  the  sound  wave,  and  R__^.  is 

2-4  expt 

the  transition  rate.  Liebermann  implies  that  his  expression  Eq.  (48), 
has  some  connection  with  a  similar  expression  for  the  absorption  of 
ultrasound  in  a  molecular  gas  found  in  a  paper  by  Richards^®.  After 
changing  the  notation  in  Richards'  paper,  his  expression  becomes. 


1 

8  =  — 
2V 


(C  -C  )  C' 
o 


2  f 

0)  f. 


C  C  f?-  +  (C  -C')  (C  -C')  U)^ 
pvlO  '•p'^v  ^ 


where  f^^  is  the  transition  rate  for  the  molecule  to  go  from  1 
is  related  to  our  transition  rate  R  in  Eq.  (41)  by 


0  and 


fj^Q  =  exp  (h  a/kg  T)  R 


while 


where  f^j^  is  the  transition  rate  for  the  molecule  to  go  from  0  ->  1, 

The  expression  of  Richards,  Eq.  (49) ,  was  first  developed  for  a  two- 

state  system  vdiere  the  relaxation  time  x  is  related  to  the  transition 
rates  in  the  usual  way,  viz.. 


W.T.  Richards,  "Supersonic  Phenomena",  Revs,  Modem  Phys ,  j^,  36 
(1939) . 
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Table  6.  Data  used  to  calculate  the  experimental  transition  rate  from  Eq.  (48)  or  Eq.  (49)  , 
usually  under  the  assianption  that  »  u. 
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(S2) 


1 

T 


10 


+  £ 


01 


However,  Richards^^  then  noted  that  must  be  neglected  when  terns 
involving  expC-2hv_,/k„T)  are  dropped,  i.e., 

Ui  D 


1 


T 


(53) 


Richards'  expression  also  agrees  with  an  expression  for  6  found  in 
Eq.  (16-15)  and  Eq.  (11-19)  in  the  book  by  Herzfeld  and  Litovitzl', 
provided  that  x  is  set  equal  to  I/^iq*  Although  other  approximations 
for  the  relaxation  time  can  be  used^in  Eq.  (49)  we  have  assumed  upon 
reading  Liebermann's  papers  that  R  ^  1^  Eq.  (48)  is 


expt 


(54) 


because  of  the  similarity  between  Eq.  (48)  and  Eq.  (49).  We  prefer 
Richards'  Eq.  (49)  over  Liebermann's  Eq.  (48)  because  a  derivation 
is  given  of  the  £ormerl6>17  whereas  the  latter  is  stated  with  no  deri- 
vation2”4^  Table  6  all  the  information  necessary  to  calculate 
an  experimental  value  of  f^^^  from  either  Eq.  (48)  or  Eq.  (49)  is  given. 

With  the  exception  of  the  specific  heats  for  benzene,  which  were 
graphically  extrapolated  from  a  paper  by  Lord,  Ahlberg,  and  Andrews  , 
all  the  data  in  Table  6  are  taken  from  the  papers  of  the  previous 
investigators.  Finally,  we  have  used  our  Eq.  (41)  to  calculate  R> 

and  thereby  f^^^  and  have  compared  it  to  the  experimental  value  of  fj^^ 

from  Liebermann's  Eq.  (48)  and  Richards'  Eq.  (49)  in  Table  7. 

CXir  results  are  in  good  agreement  with  experiment , considering  the 
approximations  made  in  our  theory  as  well  as  in  the  experiments.  Inspec¬ 
tion  of  Table  7  shows  better  agreement  of  our  theoretical  transition 
rate  f^^^  with  Richards'  experimental  f^^  for  benzene  and  p-dichloro- 

benzene  than  with  Liebermann's  experimental  transition  rate.  For 
cyclohexane  the  experimental  transition  rate  of  Liebermann  and  Richards 
are  approximately  the  same,  while  Liebermann's  value  of  f^^  for  carbon 

tetrachloride  is  more  in  agreement  with  our  theory  than  Richards ' 
value.  More  investigation  is  needed  to  determine  why  Liebermann  and 
Richards  have  different  expressions  for  the  experimental  transition 
rate. 


1 7 

K.F.  Herzfeld,  and  T.A.  Litovitz,  "Absorption  and  Dispersion  of 
Ultrasonic  Waves",  (Academic  Press,  NY,  1959),  p.  77  and  p.  63 
respectively. 

1 8 

R.C.  Lord,  Jr.,  J.E.  Ahlberg,  and  D.H.  Andrews,  "Calculation  of  the 
Heat  Capacity  Curves  of  Crystalline  Benzene  and  Benzene  - 
J.  Chem.  Phys.  5,  649  (1937). 
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IV,  DISCUSSION 


We  have  assumed  a  coupling  between  three  phonons  of  lattice 
energy  and  the  0  1  transition  in  the  lowest  vibrational  mode  of  a 

molecule,  situated  in  a  simple  cubic  lattice  described  by  the  MPRN 
model.  Then,  with  the  assumption  that  the  density  of  states  can  be 
approximated  by  sl  Debye-like  distribution,  we  calculated  the  transition 

rate  R  for  a  molecule  to  absorb  energy  from  the  lattice.  We  made 
several  approximations  in  order  to  calculate  the  transition  rate  for 
benzene,  p-dichlorobenzene,  carbon  tetrachloride,  and  cyclohexane. 

The  most  notable  are  that  the  actual  crystal  structures  can  be 
approximated  by  a  simple  cubic  lattice,  and  that  the  non-central  force 
constant  is  determined  artificially  by  making  it  consistent  with  the 
known  Debye  frequency.  We  compared  our  f^^^  against  the  experimental 

value  f^^  of  the  other  investigators,  which  was  determined  from  an 
10  2  ^ 
expression  for  absorption  in  gases  by  Liebermann  and  Richards 
This  expression  does  not  account  for  the  differences  in  absorption 
actually  measured  along  the  a,b,c  principal  crystal  axes  by  Victor 
and  Beyer^.  A  better  expression  for  absorption  must  be  developed. 

When  we  tried  to  compare  our  theoretical  transition  rate  with 
the  theoretical  one  developed  by  Liebermann^ we  realized  that  his 
final  expression  was  incorrect  because  it  failed  to  conserve  energy. 
Under  these  circumstances  we  decided  that  a  comparison  of  the  two 
expressions  was  inappropriate,  even  though  the  theoretical  transition 
rate  developed  by  Liebermann  gives  a  reasonable  value  for  benzene  of 

Our  paper  suggests  several  directions  for  future  research.  V/ith 
the  expressions  that  we  have  developed,  a  theoretical  curve  of 
absorption  versus  ultrasonic  frequency  or  temperature  can  be  plotted 
and  compared  with  the  corresponding  experimental  curve.  Many  experi¬ 
mentalists  have  suggested  that  the  exact  expression  for  the  density  of 
states  g(w.(k))  might  be  important,  instead  of  the  Debye  approximation 
Since  the  Jxact  density  of  states  is  known  for  the  MPRN^  ^12  model, 
the  differences  in  the  theoretical  transition  rate  can  be  compared  for 
this  model  using  our  formulas.  Finally,  the  actual  crystal  structure 
should  be  put  into  the  calculation  as  well  as  the  dependence  of  the 
transition  rate  on  molecular  arrangements  along  the  different  crystal 
directions , 
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APPENDIX  A 


Evaluation  of  the  Integral  I (a, a) 
The  integral  I((X,a)  is  defined  in  Eq,  C42)  as 


/ 

O  •'tr- 


x-ci-2  y=“X+C^t-l) 


g  xyC^-x-yjdy  dx _ _ 


(A.l) 


The  Taylor  series  expansions  of  etc.,  are  given  by 


ax  - 
e  -1 


n  R  Cax)"“ 
n  Qn)!  ’ 


CA.2a) 


e“y-l) 


=  1  -  f  -  ^  (-«"  B„  ^  . 


(A.  2b) 


aCa-x-y] 

aCo-x-y) 


=  1  - 


W 

>-x-y)  _  g  (aCa-x-y)) 

2  “  Z-.  n  (2n) ! 


(A. 2c) 


where  the  B  's  are  the  Bernoulli  numbers^^,  (B,  =  1/6,  B  =  1/30, 
n  ^  ^  . 

B^  =  1/42,  B^  =  1/30,  etc).  When  these  expansions  are  substituted 

into  Eq.  CA.l),  we  can  generate  an  approximation  to  I (a, a)  in  powers 
of  a.  To  obtain  our  final  form  we  have  made  use  of  the  following 
relations 


1  1 
f  y^id-yf  dy^  f  y"(d-y)“  dy  , 

•/  J  1  J  1 


(A.  3) 


^^M.R.  Spiegel,  in  "Mathematical  Handbook  of  Formulas  and  Tables"  (Schaum's 
Outline  Series,  McGraw-Hill,  NY,  1968),  p.  114. 
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where  d  =  a  -  x,  and 


/  /■  - 

•'x  =  a-2  *'y=d-l 


y"*Cd-y)^<iy<ix  = 


•/x=a.-2  Jy=d-1 


CA.4) 

Cd-y)^dydx 


Rewriting  1(0., a),  with  the  help  of  Eq.  (44)  we  have 
1  1  “ 


(A.  5) 


ICa,a)  =  /*  /*  2  (£^(x,y)a’')dydx  = 

*4=0-2  ‘4=d-l  n=0  n=0 


where  the  integrand  fj^(x,y)  obeys  the  following  relation 


(3-0)4 


■  ■/ ,  I 


fjjCx,y)dydx 


(A,  6) 


x=o-2  ‘'y=d-l 


We  list  the  f^’s  in  Table  I.  The  c^'s  have  already  been  given  as  a 
function  of  o  in  Table  1. 
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Table  I.  The  integrand  £  (x,y)  corresponding  to  the  constant  c  for  the  cases 
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APPENDIX  B 


Taylor  Series  Expansion  of 
About  the  point  a  =  3. 


The  integral  1(0-, o)  is  defined  in  Eq.  (45)  as 


(B.l) 


I (a, a)  *  a^(3-a)^ 


1  1 

f  ufp^•fl-p)u)  /  ((l-p)uv-(l-p)u-^l) 

1.0  Jy.0  (.“‘'*■’’>“'"'‘■‘”“*'’-1) 


dvdu  t 


vdiere 


p  =  a  -  2 


CB.2) 


(B.3) 


3,  is 


0<p<l>  i.e.,  2<a.<3  . 

The  Taylor  series  expansion  of  I(<x,o)  about  the  point  n 


e^°‘f2a^  -  4a  -I-  1)  *  e”Ca^  *  4a  -  2)  l]  _  ^^2 

8Ce«-  1)' 

Ff^^ri  -1  -3-30a-15a^-18a^) ■»-e”C3-lSa-15a^-3a  )  -l][a-3)  ^ 


120Ce“  -  1)^ 


+  .  •  • 


(B.4) 


where  we  have  considered  a  to  be  fixed. 
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1  University  of  California 
San  Diego 
Revelle  College 
Department  of  Physics 
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